* set path

	global path "~/agriculture_temperature"
	
	use $path/data/merged_for_reg, clear

	xtset ctycode year

* merge in GS-based bins

	merge 1:1 ctycode year using $path/data/weather_idw_1985_2015_GS
	drop _merge
	
	egen GSdd5binPos = rowtotal(GSdd5bin*)
	gen GSdd5binNeg = 245 -  GSdd5binPos
	replace GSdd5binNeg = 0 if GSdd5binNeg < 0
	drop GSdd5binPos	
	
	* ref bin
	drop GSdd5bin10
	
* results based on GS temp

	estimates clear 
	
	tempfile r1 r2 r3

********************************Value Estimates(Figure B1)***********************

	* all areas
	foreach t in 5 {
		foreach x in farming_output {
			qui reghdfe ln_`x' GSdd`t'bin* ///
					c.GSpre##c.GSpre c.GSssd##c.GSssd c.GSrhu##c.GSrhu c.GSprs##c.GSprs c.GSwin##c.GSwin ///
					i.pro#c.year##c.year ///
					if inrange(year,1995,2015) & inlist(region3,1,2), ///
					absorb(ctycode year) cluster(ctycode proXyear) 
			estimates store all_`x'_`t'C
		parmest , label saving(`"`r1'"',replace) format(p %8.2f) stars(0.1 0.05 0.01)  idstr("Aggregate")
		
	*****north*****
			qui reghdfe ln_`x' GSdd`t'bin* ///
					c.GSpre##c.GSpre c.GSssd##c.GSssd c.GSrhu##c.GSrhu c.GSprs##c.GSprs c.GSwin##c.GSwin ///
					i.pro#c.year##c.year ///
					if inrange(year,1995,2015) & inlist(region3,1), ///
					absorb(ctycode year) cluster(ctycode proXyear) 
			estimates store north_`x'_`t'C
		parmest , label saving(`"`r2'"',replace) format(p %8.2f) stars(0.1 0.05 0.01)  idstr("North")
	******south*****	
			qui reghdfe ln_`x' GSdd`t'bin* ///
					c.GSpre##c.GSpre c.GSssd##c.GSssd c.GSrhu##c.GSrhu c.GSprs##c.GSprs c.GSwin##c.GSwin ///
					i.pro#c.year##c.year ///
					if inrange(year,1995,2015) & inlist(region3,2), ///
					absorb(ctycode year) cluster(ctycode proXyear) 
			estimates store south_`x'_`t'C
		parmest , label saving(`"`r3'"',replace) format(p %8.2f) stars(0.1 0.05 0.01)  idstr("South")	

			preserve 
			use `"`r1'"',clear
			append using `"`r2'"' `"`r3'"' 
			gen crop = "`x'"
			rename id  region
			save $temp/est_cropvalGS_`x',replace
			restore			
			}
		}
		
	esttab *output*5C, b(4) se(4) keep(GSdd*bin*) brackets star(* 0.1 ** 0.05 *** 0.01)			
	
* results based on GS temp
**********************Product Estimates(Figure B2)*********************************
	tempfile r2 r3
	
	foreach x in maize  wheat rice vegetable soybean   {
		 qui reghdfe ln_`x'_yield GSdd5bin* ///
					c.GSpre##c.GSpre c.GSssd##c.GSssd c.GSrhu##c.GSrhu c.GSprs##c.GSprs c.GSwin##c.GSwin ///
			i.pro#c.year##c.year ///
			if inrange(year,1995,2015) & region3==1,  ///
			absorb(ctycode year) cluster(ctycode proXyear) 
			estimates store north_`x'_yld_5C
		parmest , label saving(`"`r2'"',replace) format(p %8.2f) stars(0.1 0.05 0.01)  idstr("North")
			
			preserve 
			use `"`r2'"',clear
			gen crop = "`x'"
			rename id  region
			save $temp/est_cropGS_north_`x',replace
			restore			
			}	
	
	*esttab north_*_yld_5C, b(4) se(4) keep(GSdd*bin*) brackets
	
	foreach x in rice vegetable oil wheat  maize  {
		 qui reghdfe ln_`x'_yield GSdd5bin* ///
					c.GSpre##c.GSpre c.GSssd##c.GSssd c.GSrhu##c.GSrhu c.GSprs##c.GSprs c.GSwin##c.GSwin ///
			i.pro#c.year##c.year ///
			if inrange(year,1995,2015) & region3==2,  ///
			absorb(ctycode year) cluster(ctycode proXyear) 
			estimates store south_`x'_yld_5C
		parmest , label saving(`"`r3'"',replace) format(p %8.2f) stars(0.1 0.05 0.01)  idstr("South")	

			preserve 
			use `"`r3'"',clear
			gen crop = "`x'"
			rename id  region
			save $temp/est_cropGS_south_`x',replace
			restore			
		}	
	
	*esttab south_*_yld_5C, b(4) se(4) keep(GSdd*bin*) brackets	

	
	use $temp/est_cropvalGS_farming_output,clear
	foreach x in maize  wheat rice vegetable soybean {
		append using $temp/est_cropGS_north_`x'
	}
	foreach x in rice vegetable oil wheat  maize {
		append using $temp/est_cropGS_south_`x'
	}	
	
	// keep temp variables only
	keep if strmatch(parm,"*dd5bin*")
	
	//level varibables
	gen level = 2.5 if strmatch(parm, "GSdd5bin0")
	replace level = 7.5 if strmatch(parm, "GSdd5bin5")
	replace level = -2.5 if strmatch(parm, "GSdd5binNeg")
	*replace level = 12.5 if strmatch(parm, "GSdd5bin10")
	replace level = 17.5 if strmatch(parm, "GSdd5bin15")
	replace level = 22.5 if strmatch(parm, "GSdd5bin20")
	replace level = 27.5 if strmatch(parm, "GSdd5bin25")
	replace level = 32.5 if strmatch(parm, "GSdd5bin30")
	replace level = 37.5 if strmatch(parm, "GSdd5bin35")

	
	// crop name 
	replace crop = "Corn" if strmatch(crop,"maize")
	replace crop = "Rice" if strmatch(crop,"rice")
	replace crop = "Wheat" if strmatch(crop,"wheat")
	replace crop = "Oil Crop" if strmatch(crop,"oil")
	replace crop = "Vegetable" if strmatch(crop,"vegetable")
	replace crop = "Soybean" if strmatch(crop,"soybean")
		
	rename estimate coef 
	rename stderr se 
	rename min95 ci_l
	rename max95 ci_h
	
	// add zero-value lines
	bysort region crop : gen last = _n == _N
	expand 2 if last
	sort region crop parm
	replace level = 12.5 if parm==parm[_n-1]
	replace coef = 0 if parm==parm[_n-1]
	replace se = 0 if parm==parm[_n-1]
	replace ci_h = 0 if parm==parm[_n-1]
	replace ci_l = 0 if parm==parm[_n-1]			
	
	save $visual/visual_croppingGS,replace	
	
	// delete dta data 
	cd $temp/
	fs *.dta
	foreach f in `r(files)'{
	   erase "`f'"
	}
	
	
